# load("./data/Robjects/Annotated_Results_LvV.RData")
# load("./data/Robjects/Ensembl_annotations.RData")
# load("./data/Robjects/DE.RData")
# load("./data/Robjects/01_DDS.RData")
# load("./data/Robjects/02_annot.RData")
load("./data/Robjects/03_DDS.RData")
for now only using the gastroc tissue
# res_sorted <- data.frame(res.gastroc, stringsAsFactors = F) %>%
# na.omit() %>%
# subset(., (log2FoldChange >= 1 |
# log2FoldChange <= -1)) %>% # set arbitrarily
# arrange(desc(abs(log2FoldChange))) # order by fold change (absolute value)
# # up-regulated genes
# up <- na.omit(res_sorted) %>%
# subset(., log2FoldChange > 0 & padj < 0.01)
res_sorted<-data.frame(gene = rownames(res.gastroc),res.gastroc,stringsAsFactors = F)%>%
na.omit()%>%
subset(., ( log2FoldChange >=1 | log2FoldChange <=-1))%>% # set abritarily
arrange(desc( abs(log2FoldChange) )) # order by fold change (absolute value)
# up-regulated genes
up<-data.frame(rn= rownames(res_sorted),res_sorted) %>%
na.omit()%>%
subset(., log2FoldChange > 0 & padj < 0.05)
gp_up <- gost(row.names(up), organism = "mmusculus", evcodes = TRUE, ordered_query = TRUE,
correction_method = "gSCS", user_threshold = 0.01,
sources = c("GO:BP", "GO:MF", "GO:CC", "KEGG", "REAC", "TF", "MIRNA", "CORUM", "HP", "HPA", "WP")) # ordered by foldchange; evcodes generates gene names (makes run time longer though)
require(enrichplot)
require(ggplot2)
require(viridis)
require(shiny)
require(DOSE)
div(gostplot(gp_up, interactive = TRUE), align = "right") # interactive plot; circle sizes are in accordance with corresponding term size; term location on the x-axis is fixed and terms from the same GO subtree are located closer to each other
# define as enrichResult object and plot
gp_up_df<-gp_up$result[,c("query", "source", "term_id", "term_name",
"p_value", "query_size", "intersection_size",
"term_size", "effective_domain_size", "intersection")]%>%
mutate(intersection = gsub(",", "/", intersection),
GeneRatio = paste0(intersection_size, "/", query_size),
BgRatio = paste0(term_size, "/", effective_domain_size))
names(gp_up_df) = c("Cluster", "Category", "ID", "Description", "p.adjust",
"query_size", "Count", "term_size", "effective_domain_size", "geneID", "GeneRatio", "BgRatio")
row.names(gp_up_df) = gp_up_df$ID
#gp_up_df<-subset(gp_up_df, Category == "REAC") # can prefilter for databases of interest
# dotplot
enrichplot::dotplot(new("enrichResult", result = gp_up_df))+
theme_minimal()+
scale_color_viridis()+
theme(axis.title=element_text(size=rel(1.5), colour = "black"),
axis.line=element_line(colour = "grey40"),
panel.grid.major=element_blank(),
panel.grid.minor=element_blank(),
axis.ticks.x=(element_line(size=0.5, colour = "grey40")))
# barplot
barplot(new("enrichResult", result = gp_up_df), showCategory = 10)+
ggplot2::facet_grid(~Cluster) +
ggplot2::ylab("Intersection size")+
theme_minimal()+
scale_fill_viridis()+
theme(axis.title=element_text(size=rel(1.5), colour = "black"),
axis.line=element_line(colour = "grey40"),
panel.grid.major=element_blank(),
panel.grid.minor=element_blank(),
axis.ticks.x=(element_line(size=0.5, colour = "grey40")))
# terms and gene names
enrichplot::heatplot(new("enrichResult", result = gp_up_df), showCategory = 10) +
theme_minimal() +
theme(
axis.title = element_text(size = rel(1.5), colour = "black"),
axis.text.x = element_text(angle = 90),
axis.line = element_line(colour = "grey40"),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.ticks.x = (element_line(size = 0.5, colour = "grey40"))
)
# terms and gene names 2
enrichplot::cnetplot(
new("enrichResult", result = gp_up_df),
showCategory = 3,
cex_label_gene = 0.5
) +
scale_colour_manual(values = c("grey70", "black"), guide = FALSE)
# terms and gene names 3
enrichplot::cnetplot(
new("enrichResult", result = gp_up_df),
showCategory = 3,
circular = TRUE,
colorEdge = TRUE,
cex_label_gene = 0.5
)